# Fit models for main and by party, without controls
main_results <- lm_robust(message_main ~ factor(frame), data)
dems_results <- lm_robust(message_main ~ factor(frame), dems)
reps_results <- lm_robust(message_main ~ factor(frame), reps)

# Collect models in a list
models <- list(main_results, dems_results, reps_results)

# Create the LaTeX table using the original models but overriding the standard errors and coefficients
tab <- texreg(models,
               custom.model.names = c("Full sample", "Democrats", "Republicans"),
               custom.coef.map = list("factor(frame)blame" = "Climate Justice",
                                      "factor(frame)global" = "Global Solidarity",
                                      "factor(frame)none" = "Appeal Only",
                                      '(Intercept)' = 'Constant'),
               caption = "Effect of climate emissions frames on Aid Support Index",
               label = "tab:aid.index",
               include.rmse = FALSE,
               include.constant = TRUE, # Ensure constant is included
               include.ci = FALSE,
               stars = c(0.01, 0.05, 0.1),
               float.pos = "h")
tab
# Write the LaTeX table to a file
writeLines(tab, "MessagingFigsR/main_reg_RR.tex")


# p-value randomization inference
# Run the randomization inference for 'none' vs 'blame'
index_none_blame <- randomization_inference(data, "message_main", "frame", "none", "blame")

index_none_global <- randomization_inference(data, "message_main", "frame", "none", "global")

sec_none_blame <- randomization_inference(data, "c_donate_1", "frame", "none", "blame")

sec_none_global <- randomization_inference(data, "c_donate_1", "frame", "none", "global")

# Create a data frame to store the results
p_value_table <- data.frame(
  Variable = c("Aid Index", "Aid Index", "Personal Donation", "Personal Donation"),
  Comparison = c("Control vs. Justice", "Control vs. Solidarity", 
                 "Control vs. Justice", "Control vs. Solidarity"),
  P_Value = c(index_none_blame$p_value, 
              index_none_global$p_value, 
              sec_none_blame$p_value, 
              sec_none_global$p_value)
)

# Rename columns as needed
colnames(p_value_table) <- c("Variable", "Comparison", "P-Value")

# Print the table (use xtable if you want LaTeX output)
print(xtable(p_value_table, include.rownames = FALSE))


# Determine unique variables across all models
unique_variables <- unique(unlist(lapply(models, function(mod) grep("^factor\\(frame\\)", names(coef(mod)), value = TRUE))))

# Extract coefficients and standard errors
results <- extract_coefficients(models, unique_variables)
results
results$model <- c('main_j', 'main_s','dem_j', 'dem_s','rep_j','rep_s')
results <- na.omit(results)
results$column_set <- factor(c('All', 'All', 'Democrats', 'Democrats', 'Republicans', 'Republicans'), levels = c("All", "Democrats", "Republicans"))

results$coefficients <- as.numeric(results$coefficients)
results$standard_errors <- as.numeric(results$standard_errors)

# Plot results
plot_results(
  results, 
  x_labels = c("main_j" = "Justice", "main_s" = "Solidarity", "dem_j" = "Justice", "dem_s" = "Solidarity", "rep_j" = "Justice", "rep_s" = "Solidarity"), 
  y_label = "Estimated Effect of Message Frame on Aid Support Index", 
  colors = c("All" = 'purple3', 'Democrats' = 'blue3', 'Republicans' = 'red3'), 
  color_column = 'column_set',
  file_path = "MessagingFigsR/main_reg_RR.png"
)

# Fit models for private donations outcome using the refactored function
sec_results <- lm_robust(c_donate_1 ~ factor(frame), data)
sec_dems_results <- lm_robust(c_donate_1 ~ factor(frame), dems)
sec_reps_results <- lm_robust(c_donate_1 ~ factor(frame), reps)

# List of models for private donations
sec_models <- list(sec_results, sec_dems_results, sec_reps_results)
secp_values <- lapply(sec_models, function(model) summary(model)$coefficients[, "Pr(>|t|)"])

# Create the LaTeX table using the original models but overriding the standard errors and coefficients
tab2 <- texreg(sec_models,
               custom.model.names = c("Full sample", "Democrats", "Republicans"),
               custom.coef.map = list("factor(frame)blame" = "Climate Justice",
                                      "factor(frame)global" = "Global Solidarity",
                                      "factor(frame)none" = "Appeal Only",
                                      '(Intercept)' = 'Constant'),
               caption = "Effect of climate emissions frames on personal donation attitudes",
               label = "tab:personal.donation",
               include.rmse = FALSE,
               include.constant = TRUE, # Ensure constant is included
               include.ci = FALSE,
               stars = c(0.01, 0.05, 0.1),
               float.pos = "h")
tab2
# Write the LaTeX table to a file
writeLines(tab2, "MessagingFigsR/sec_reg_RR.tex")

# Determine unique variables across all private donations models
sec_unique_variables <- unique(unlist(lapply(sec_models, function(mod) grep("^factor\\(frame\\)", names(coef(mod)), value = TRUE))))

# Extract coefficients and standard errors for private donations models
sec_results_df <- extract_coefficients(sec_models, sec_unique_variables)
sec_results_df
sec_results_df$model <- c('main_j', 'main_s', 'dem_j', 'dem_s', 'rep_j', 'rep_s')
sec_results_df <- na.omit(sec_results_df)
sec_results_df$column_set <- factor(c('All', 'All', 'Democrats', 'Democrats', 'Republicans', 'Republicans'), levels = c("All", "Democrats", "Republicans"))

sec_results_df$coefficients <- as.numeric(sec_results_df$coefficients)
sec_results_df$standard_errors <- as.numeric(sec_results_df$standard_errors)

# Plot results for private donations
plot_results(
  sec_results_df, 
  x_labels = c("main_j" = "Justice", "main_s" = "Solidarity", "dem_j" = "Justice", "dem_s" = "Solidarity", "rep_j" = "Justice", "rep_s" = "Solidarity"), 
  y_label = "Estimated Effect of Message Frame on Personal Donations", 
  colors = c("All" = 'purple3', 'Democrats' = 'blue3', 'Republicans' = 'red3'), 
  color_column = 'column_set',
  file_path = "MessagingFigsR/sec_reg_RR.png"
)

# Define the variables for mechanisms
outcomes <- c("collective_guilt", "social_prox", "worriedwarming")
parties <- c('Democrat', 'Republican')

# Create an empty list to store regression results for mechanisms
reg_results <- list()

for (party in parties) {
  party_data <- data[data$party_final_analysis == party, ]
  for (outcome in outcomes) {
    reg_model <- lm_robust(as.formula(paste(outcome, "~ factor(frame)")), party_data)
    reg_results[[paste(outcome, party, sep = "_")]] <- reg_model
  }
}

# Fit models for the entire sample for each outcome
collective_guilt_All <- lm_robust(collective_guilt ~ factor(frame), data)
social_prox_All <- lm_robust(social_prox ~ factor(frame), data)
worriedwarming_All <- lm_robust(worriedwarming ~ factor(frame), data)

# List of mechanism models
mech_models <- list(collective_guilt_All, reg_results$collective_guilt_Democrat, reg_results$collective_guilt_Republican, 
                    social_prox_All, reg_results$social_prox_Democrat, reg_results$social_prox_Republican,
                    worriedwarming_All, reg_results$worriedwarming_Democrat, reg_results$worriedwarming_Republican)

mech_models2 <- list(collective_guilt_All, reg_results$collective_guilt_Democrat, reg_results$collective_guilt_Republican, 
                    social_prox_All, reg_results$social_prox_Democrat, reg_results$social_prox_Republican)

mech_models3 <- list(worriedwarming_All, reg_results$worriedwarming_Democrat, reg_results$worriedwarming_Republican)

# Create the LaTeX table using the original models but overriding the standard errors and coefficients
tab3 <- texreg(mech_models2,
               custom.model.names = c("Guilt Full", "Guilt D", "Guilt R", 'Proximity Full','Proximity D','Proximity R'),
               custom.coef.map = list("factor(frame)blame" = "Climate Justice",
                                      "factor(frame)global" = "Global Solidarity",
                                      "factor(frame)none" = "Appeal Only",
                                      '(Intercept)' = 'Constant'),
               caption = "Effect of message frame on collective guilt and social proximity",
               label = "tab:mechanisms",
               include.rmse = FALSE,
               include.constant = TRUE, # Ensure constant is included
               include.ci = FALSE,
               stars = c(0.01, 0.05, 0.1),
               float.pos = "h")
tab3
# Write the LaTeX table to a file
writeLines(tab3, "MessagingFigsR/mechanisms_RR.tex")

tab4 <- texreg(mech_models3,
               custom.model.names = c("Worried All", "Worried D", "Worried R"),
               custom.coef.map = list("factor(frame)blame" = "Climate Justice",
                                      "factor(frame)global" = "Global Solidarity",
                                      "factor(frame)none" = "Appeal Only",
                                      '(Intercept)' = 'Constant'),
               caption = "Effect of message frame on concern about global warming",
               label = "tab:mechanisms_warming",
               include.rmse = FALSE,
               include.constant = TRUE, # Ensure constant is included
               include.ci = FALSE,
               stars = c(0.01, 0.05, 0.1),
               float.pos = "h")
tab4
# Write the LaTeX table to a file
writeLines(tab3, "MessagingFigsR/app_mech_worry_RR.tex")


# Determine unique variables across all mechanism models
mech_unique_variables <- unique(unlist(lapply(mech_models, function(mod) grep("^factor\\(frame\\)", names(coef(mod)), value = TRUE))))

# Extract coefficients and standard errors for mechanism models using HC2 standard errors
mech_results_df <- extract_coefficients(mech_models, mech_unique_variables)
mech_results_df$model <- c('guilt a j', 'guilt a s', 'guilt d j','guilt d s','guilt r j','guilt r s',
                          'prox a j','prox a s','prox d j','prox d s','prox r j','prox r s',
                           'worry a j','worry a s','worry d j','worry d s','worry r j','worry r s')

# Clean up the results data frame for mechanisms
mech_results_df <- na.omit(mech_results_df)
mech_results_df$column_set <- factor(c('Collective Guilt- Justice','Collective Guilt- Justice','Collective Guilt- Justice','Collective Guilt- Justice','Collective Guilt- Justice','Collective Guilt- Justice',
                                        'Social Proximity- Solidarity','Social Proximity- Solidarity','Social Proximity- Solidarity','Social Proximity- Solidarity','Social Proximity- Solidarity','Social Proximity- Solidarity',
                                       'Concern over Global Warming','Concern over Global Warming','Concern over Global Warming','Concern over Global Warming','Concern over Global Warming','Concern over Global Warming'), 
                                     levels = c('Collective Guilt- Justice','Social Proximity- Solidarity','Concern over Global Warming'))
mech_results_df$party <- c('All','All','Democrats','Democrats','Republicans','Republicans',
                           'All','All','Democrats','Democrats','Republicans','Republicans',
                           'All','All','Democrats','Democrats','Republicans','Republicans')

mech_results_df$coefficients <- as.numeric(mech_results_df$coefficients)
mech_results_df$standard_errors <- as.numeric(mech_results_df$standard_errors)

# Remove unnecessary models and columns
mech_results_df <- mech_results_df[!(mech_results_df$model %in% c('guilt a s','guilt d s','guilt r s','prox a j','prox d j','prox r j')), ]
mech_results_df <- mech_results_df[!(mech_results_df$column_set == 'Concern over Global Warming'),]

plot_results(
  mech_results_df, 
  x_labels = c("guilt d j" = "", "guilt a j" = "", "guilt d s" = "Solidarity", "guilt r j" = "", "guilt r s" = "Solidarity",
               "prox d j" = "Justice", "prox d s" = "", "prox a s" = "", "prox r j" = "Justice", "prox r s" = "",
               "worry a j" = "Justice", "worry a s" = "Solidarity", "worry d j" = "Justice", "worry d s" = "Solidarity",
               "worry r j" = "Justice", "worry r s" = "Solidarity"), 
  y_label = "Estimated Effect of Message Frame on Attitudes", 
  color_column = 'party',
  colors = c('All'= 'purple','Democrats' = 'blue3', 'Republicans' = 'red3'), 
  file_path = "MessagingFigsR/mechanisms_RR.png"
)

# Compute results table for t-tests
control_group <- subset(data, factor(frame) == 'none')
control_group$group <- ifelse(control_group$white == 1, "White",
                              ifelse(control_group$black == 1, "Black", NA))

t_test_message_main <- t.test(message_main ~ group, data = control_group, subset = !is.na(group))
t_test_c_donate_1 <- t.test(c_donate_1 ~ group, data = control_group, subset = !is.na(group))

results_table <- create_results_table(t_test_message_main, t_test_c_donate_1, control_group)
colnames(results_table) <- c("Outcome", "Mean Difference", "t Statistic", "p-value", "95% CI Lower", "95% CI Upper", "Mean Black", "Mean White")
latex_table <- xtable(results_table, caption = "Generosity Comparison Control Blacks vs. Whites", label = "tab:control.generosity")

# Alternatively, save the table to a .tex file
writeLines(print(latex_table, include.rownames = FALSE), "MessagingFigsR/control_generosity_RR.tex")


